% Read Agilent fid file, and then rip off only the data block of each
% channel, then save the data to variables respectively.
% The output has real and imaginary number both.

function [realData, imaginaryData, headerInfo] = readfid(filename, channel)

%     s_fileHeader=32;  % the bytes of file header %
%     s_blockHeader=28; % the bytes of block header %
   

    % Read fid file and save the binary information %
    fileID = fopen(filename, 'r', 'b');
    
    %% The header information %%
    % the file header information %
    nblocks   = fread(fileID,1,'int32');
    ntraces   = fread(fileID,1,'int32');
    np        = fread(fileID,1,'int32');
    ebytes    = fread(fileID,1,'int32');
    tbytes    = fread(fileID,1,'int32');
    bbytes    = fread(fileID,1,'int32');
    vers_id   = fread(fileID,1,'int16');
    status    = fread(fileID,1,'int16');
    nbheaders = fread(fileID,1,'int32');
    dataPoints = np/2;
        
    % status information %
    s_data    = bitget(status,1);
    s_spec    = bitget(status,2);
    s_32      = bitget(status,3);
    s_float   = bitget(status,4);
    s_complex = bitget(status,5);
    
    % make empty array for save data %
    realData=zeros(dataPoints, 1);
    imaginaryData=zeros(dataPoints, 1);
    
    % error checker %
    if channel > nblocks 
        fprintf('[Error] Incorrect channel number : this data file has only %d channel \n', nblocks);
        return;
    end
    
    % calculating the seeking position %
    seek_indx = (channel-1)*bbytes;
        
    % Seeking the file for reading data %
    s_fseek = fseek(fileID, seek_indx, 'cof');
    
    % error check of the fseek operation %
    if s_fseek == -1
        return;
    end    
       
    % the block header information %
    scale     = fread(fileID,1,'int16');
    bstatus   = fread(fileID,1,'int16');
    index     = fread(fileID,1,'int16');
    mode      = fread(fileID,1,'int16');
    ctcount   = fread(fileID,1,'int32');
    lpval     = fread(fileID,1,'float32');
    rpval     = fread(fileID,1,'float32');
    lvl       = fread(fileID,1,'float32');
    tlt       = fread(fileID,1,'float32');

    % get complex data %
    for idx2 = 1 : dataPoints
        realData(idx2) = fread(fileID, 1, 'float32');
        imaginaryData(idx2) = fread(fileID, 1, 'float32');    
    end        
    
    % Header information of file and blocks %
    headerInfo = [nblocks, ntraces, np, ebytes, tbytes, bbytes, vers_id, status, nbheaders];
    headerInfo = [headerInfo, scale, bstatus, index, mode, ctcount, lpval, rpval, lvl, tlt];

    %% closing file %%   
    fclose(fileID);    
    return;

end  % the last line of function %